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Drag reduction in stationary turbulent fiows by bubbles is sensitive to the dynamics of bubble 
oscillations. Without this dynamical effect the bubbles only renormalize the fluid density and 
viscosity, an effect that by itself can only lead to a small percentage of drag reduction. We show in 
this paper that the dynamics of bubbles and their effect on the compressibility of the mixture can 
lead to a much higher drag reduction. 
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I. INTRODUCTION 

Drag reduction in turbulent flows is a subject of tech- 
nological importance and of significant basic interest. As 
is well known, drag reduction can be achieved using a 
number of additives, including flexible polymers, rod-like 
polymers and fibers, surfactants, and bubbles [1]. While 
the subject of drag reduction by polymers had seen rapid 
theoretical progress in the last few years [2-7] the under- 
standing of drag reduction by bubbles lags behind. For 
practical applications in the shipping industry the use of 
polymers is out of the question for economic and envi- 
ronmental reasons, but air bubbles are potentially very 
attractive. 

The theory of drag reduction by small concentrations 
of minute bubbles is relatively straightforward, since un- 
der such conditions the bubbles only renormalize the 
density and the viscosity of the fluid, and a one-fluid 
model suffices to describe the dynamics [8]. The fluid 
remains incompressible, and the equations of motion are 
basically the same as for a Newtonian fluid with renor- 
malized properties. The amount of drag reduction under 
such conditions is however limited. But when the bub- 
bles increase in size, the one-fluid model loses its validity 
since the bubbles become dynamical in the sense that 
they are no longer Lagrangian particles, their velocity 
is no longer the fluid velocity at their center, and they 
begin to fluctuate under the influence of local pressure 
variations. The fluctuations of the bubbles are of two 
types: 1) the bubbles are no longer spherical, distorting 
their shape according to the pressure variations, and 2) 
the bubbles can oscillate radially (keeping their spherical 
shape) due to the commpressiblity of the gas inside the 
bubble. The first effect was studied numerically using the 
"front tracking" algorithm in Ref [9, 10]. However, the 
results indicate either a drag enhancement, or a limited 
and transient drag reduction. This leads one to study the 
possibility of explaining bubbly drag reduction by bub- 
ble oscillations. Indeed, a theoretical model proposed by 
Legner [11] successfully explained the bubbly drag reduc- 
tion by modifying the turbulent viscosity in the bubbly 
flow by the bulk viscosity of the bubbles. While the bulk 
viscosity is important only when the bubbles are com- 
pressible, it is important and interesting to see how and 
why it affects the charactistics of the flow. The aim of 
this paper is to study the drag reduction by bubbles when 
bubble oscillations are dominant. Finally we compare our 



finding with the results in Ref [11], showing that a non- 
physical aspect of that theory is removed, while a good 
agreement with experiment is retained. 

In our thinking we were influenced by two main find- 
ings, one experimental and the other simulational. The 
experiment [12] established the importance of bubble dy- 
namics in effecting drag reduction. The same turbulent 
flow was set up once in the presence of bubbles and once 
in the presence of glass spheres whose density was smaller 
than that of the ambient fluid. While bubbles effected 
drag reduction for sufficiently high Reynolds number, the 
glass spheres enhanced the drag. In the simulation [13] it 
was demonstrated that the drag reduction by the bubbles 
is connected in an intimate way to the effective compress- 
ibility of the mixture. (The fiuid by itself was taken as 
incompressible in the simulation). These two observa- 
tions, in addition to the experiments [14] will be at the 
back of our mind in developing the theory, with the final 
elucidation of all these observations in the last sections 
of this paper. 

Im Sect. II we present the average (field) equations for 
fluids laden with bubbles. This theory follows verbatim 
earlier work [15-17] and it is limited to rather small bub- 
bles (of the order of the Kolmogorov scale) and to poten- 
tial flows. In Sect. Ill we employ the theory to find out at 
which Reynolds and Weber numbers the bubbles interact 
sufficiently strongly with the fluid to change significantly 
the stress tensor beyond simple viscosity renormalization. 
In Sect. IV we study the balance equation for momen- 
tum and energy in the turbulent boundary layer. This 
leads to the main section of this paper. Sect. V which 
presents the predictions of the theory regarding drag re- 
duction by bubbles. The oscillations of the bubbles at 
sufficiently high Weber numbers are shown to be an im- 
portant physical reason for the phenomenon. A summary 
and discussion are presented in Sect. VI. 



II. AVERAGED EQUATIONS FOR BUBBLY 
FLOWS 



A Newtonian fluid with density p is laden with bubbles 
of density , and radius R which is much smaller than 
the outer scale of turbulence C. The volume fraction of 
bubbles C is taken sufficiently small such that the direct 
interactions between bubbles can be neglected. In writ- 
ing the governing equations for the bubbly fiow we will 
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assume that the length scales of interest are larger than 
the bubble radius. Later we will distinguish however be- 
tween the case of microbubbles whose radius is smaller 
than the Kolmogorov scale t] and bubbles whose radius 
is of the order of rj or slightly larger. For length scales 
larger than the bubbles one writes [15-17]: 



equation of motion for each bubble 



-CVp + CV -a-CF, 



(1) 



9/i 



(U-w) 



p fDU 
2\Dt 



where fi is the dynamical viscosity of the neat fluid. 
In this equation the force acting on the bubble is 
only approximate, since we neglect gravity, the lift 
force and the add-mass force due bubble oscilla- 
tions. We include only the viscous force and the 
add-mass force due to bubble acceleration, and we 
will show that this is sufficient for enhancing the 
drag reduction by the bubble dynamics. It can be 
argued that adding the other forces does not change 
things qualitatively. 

• equation of motion for the carrier fluid 

p{l-C)^ = -(1-C) Vp+(1-C) V-a+CF+CV-T ; 

(2) 



continuity equation 
d{l-C) 



dt 



+ V • (1 - C) J7 = . 



(3) 



In these equations, U and w is the velocity of the carrier 
fluid and the bubble respectively, and 

F and r are the force and the stress caused by the dis- 
turbance of the flow due to the bubbles, the Lagrangian 
derivatives are defined by 



and 



Da da 



da „ 

a = — + w ■ Va 
ot 



(5) 



(6) 



As the density of the bubble is usually much smaller than 
the fluid, is taken to be 0. Combining (1) and (2), we 
have 



p{l - C) 



DU 



-Vp+V -a + CV -T . (7) 



Note that the term containing F disappears in the last 
equation because of the cancellation of action and reac- 
tion forces. 

The bubbles affect the flow in two ways: 



• changing the effective density of the fluid; 

• introducing an additional stress tensor r to the 
fluid velocity equation (7). 

The expression used for r is extremely important for the 
discussion at hand. It is commonly accepted that the 
stress tensor is affected by three factors: 



(8) 



In this equation is the viscous stress tensor, written 
as: 



5m (dlh 
2 \dxi 



dxi 



(9) 



For very small bubbles (micro-bubbles) of very small den- 
sity this is the only significant contribution in Eq. (8). 
When this is the case; the bubble contribution to the 
stress tensor can be combined with a in Eq. (7), resulting 
in the effective viscosity given by 



Meff = (1 + -C)m • 



(10) 



The study of drag reduction under this renormalization 
of the viscosity and the density was presented in Ref. [8] , 
with the result that drag reduction can be obtained by 
putting the bubbles out of the viscous sub-layer and not 
too far from the wall. The amount of drag reduction is 
however rather limited in such circumstances. 

The other two contributions in Eq. (8) arc the concern 
of the present paper. The component is non zero only 
when the bubble is not a Lagrangian particle, having a 
relative velocity w U with respect to the fluid; then 
the bubble radius is changing in time. Explicitly [15-17]: 



R^ + ^{w-U)-{w-U) 



in) 



The last contribution Tg is sensitive to the change in pres- 
sure of the fluid due to the bubbles. It reads [15]: 

T,=~ J{p-p°)nndA. (12) 

Here p^ is the pressure of the fluid without bubbles, n 

is the normal unit vector to the bubble surface, and dA 
is the area differential. The relation of this expression to 
the relative velocity and to the bubble dynamics calls for 
a calculation, which in general is rather difficult. Such 
a calculation was achieved explicitly only for potential 
flows, with the final result [15, 16]: 



"^{w-u) ■{w-u)-RR~ hkf 
-g(«, -[/)(«,-[/). 



(13) 
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III. RELATIVE IMPORTANCE OF THE 
STRESS CONTRIBUTIONS AS A FUNCTION OF 
THE REYNOLDS NUMBER 

The relative importance of the three contributions Tiy, 
and Tg depends on the Reynolds number and on R/jC. 
To study this question represent Eq. (1) as follows 



with 



DU 



w 



-Vp + V • a- 



9fi 



w) 



(14) 



Consider first the case of small bubble size, R < r], and 
small Reynolds numbers. In this case the viscous term 
on the RHS is dominant, and the difference between U 
and w cannot be large. The bubbles behave essentially 
as Lagrangian tracers. On the other hand, at high values 
of Re and for larger bubbles, R> rj, the term Vp should 
be re-interpreted on the scale of the bubble as 



Vp 



= P- 



p{x + R)- p{x - R) 
2R 



(15) 



2R 



where x is the location of the bubble. The second line in 
Eq. (15) follows from Bernouli's equation. When the size 
of the bubble becomes of the order of the Kolmogorov 

scale or larger, wc have 



U^{x + R)-U^{x-R) ~ 2U{x){eR) 

~ 2U{x)U,ms 



1/3 



(16) 



1/3 



where C/rms is the r.m.s. of the turbulent velocity. At 
this point we can ask what is the value of the Reynolds 
number for which the viscous term is no longer domi- 
nant, allowing for significant fiuctuations inU ~w. This 
happens when the terms in Eq. (14) are comparable, i.e. 
when 



\U-w\ ~ 



R 

C/^sRe (R 



1/3 



i?2 



(17) 



4/3 



This equation contains an important prediction for ex- 
periments. It means that the fluctuations in the relative 
velocity of the bubble with respect to fluid is of the or- 
der of the outer fluid velocity when Re is larger than 
[R/C)"^/^. In most experiments, R/C - 0{lQ-'^) and it 
is therefore sufficient to reach Re ~ O(IO^) for |f/ — tu| 
to be of the order of ?7rms- Note that this is precisely the 
result of the experiment [12]. 

This discussion has consequences for the bubble dy- 
namics and oscillations. At small Re, w — U\s small and 
T Tv Then the equation of the mixture becomes: 



Peff 



DU 

W 



(18) 



Peff = p(l-C), 

5 

<Teff = <t(1 + - C) , 



(19) 

(20) 



meaning that only the effective density and viscosity are 
changed, as is usually assumed in numerical simulations 
of "point" bubbles [12, 13]. On the other hand, when Re 
is large \w — U\ is comparable to J7rms- This will affect 
the stress tensor on scales larger than the bubble size via 



Tj^ and Tg. Furthermore, 



RR+-ie 



27 
R 



-{w-U)-{w-U)~4:ii 



R 
R 
(21) 

where 7 is the surface tension. This equation tells us that 
the radial oscillations of the bubbles are excited by the 
relative velocity w — U. When w — U = 0, then 



27^ 



(22) 



and so i? is a constant. Similarly, R is small ii w — U 
is small. The strength of the oscillation can be charac- 
terised by the Weber number 



We = 



p\w-U\'^R 
7 



(23) 



As a summary, the additional stress tensor r in the 
basic Eq. (7) due to the presence of bubble is a sum of 
three contributions, Ti,,t^, and r^, see Eq. (8). By using 
Eqs. (9), (11) and (13), we have 



= p{ [i {w-U)-iw- U) 



1 
4 



RR 



R' 



(24) 



where the tensor S has only one nonzero component, 
Sxy = S. The relative importance of the various terms 
in T depends on the values of Re and We. If We is suffi- 
ciently large, there will be a large change in the diagonal 
part of Tg . In the following section we show that this can 
be crucial for drag reduction. 



IV. BALANCE EQUATIONS IN THE 
TURBULENT BOUNDARY LAYER 

At this point we apply the formalism detailed above 
to the question of drag reduction by bubbles in a sta- 
tionary turbulent boundary layer with plain geometry. 
This can be a pressure driven turbulent channel flow or a 
plain Couctte flow, which is close to the circular Couctte 
flow realized in [12]. Let the smallest geometric scale 
be 2L (for example the channel height in a channel flow) , 
the unit vector in the streamwise and spanwise directions 
be X and z respectively, and the distance to the nearest 
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wall he y L. The velocity U{r, t) has only one mean 
component, denoted hy V = xV, that depends only on 
y: V = V{y). Denoting turbulent velocity fluctuations 
(with zero mean) by u(r, t) we have the Reynolds decom- 
position of the velocity field to its mean and fluctuating 
part: 

C/(r,t) = y(y)x + w(r,i) . (25) 

Long time averages are denoted by (. . . ). Having dynam- 
ical equations (7) and (24), we can consider the effect of 
the bubbles on the statistics of turbulent channel flow. 
For this goal we shall use a simple stress model of planar 
turbulent flow. A similar model was successfully used in 
the context of drag reduction by polymeric additives [18]. 
This model is based on the balance equations of mechan- 
ical momentum, which we consider in the next Sec. IV A 
and the balance of the turbulent kinetic energy, discussed 
in Sec. IV B. The variables that enter the model are the 
mean shear 

S = dV/dy , (26a) 
the turbulent kinetic energy density 

^=^P(1-C)(H2), (26b) 

and the Reynolds stress 

W^y = -p{l - C){u^Uy) . (26c) 

A. Momentum balance 

From Eq. (7) we derive the exact equation for the 
momentum balance by averaging and integrating in the 
usual way, and find fov y <^ L: 

P = flS + W^y + C{T^y) . (27) 

Here P is the momentum flux toward the wall. In a 

channel flow P = p'L, were p' = —dp/dx is the (con- 
stant) mean pressure gradient. In a plain Couette flow P 
is another constant which is determined by the velocity 
difference between the two walls. For C = Eq. (27) is 
the usual equation satisfied by Newtonian fluids. 

To expose the consequences of the bubbles we notice 
that the diagonal part of the bubble stress tensor r [the 
first line in the RHS of Eq. (24)] does not contribute to 
Eq. (27). The xy component of the off-diagonal part 
of T is given by the 2nd line in Eq. (24). Define the 
dimensionless ratio 

^ = {{Wx - Ux){Wy - Uy)) 
2{UxUy) 

For later purposes it is important to assess the size and 
sign of a. For small values of Re, a. is small according 
to Eq. (29). On the other hand, it was argued in [19, 



20] that for large Re the fluctuating part of w is closely 
related to the fluctuating part of u. The relation is 

w-U ^2u . (29) 

If we accept this argument verbatim this would imply 
that a w 2 and is positive definite, as we indeed assume 
bellow. With this definition we can simplify the appear- 
ance of Eq. (27): 

P = Meff^ + (30) 

with /ieir defined by Eq. (10). Below we consider the high 
Re limit, and accordingly can neglect the first term on 
the RHS. 



B. Energy balance 

Next, we consider the balance of turbulent energy in 
the log-layer. In this region, the production and dissipa- 
tion of turbulent kinetic energy is almost balanced. The 
production can be calculated exactly, WxyS. The dissi- 
pation of the turbulent energy is modelled by the energy 
fiux which is the kinetic energy K{y) divided by the typ- 
ical eddy t urn over t ime at a distance y from the wall, 
which is ■y/p(l — C)y/b\/K where 6 is a dimensionless 
number of t he order o f unity. Thus the fiux is written 
as hK'^l'^ jy^ p(\ — C). The extra dissipation due to the 
bubble is C{TijSij) where Sij = dui/dxj. In summary, 
the turbulent energy balance equation is then written as: 

I + C{Tii Sij) = W^yS . (31) 

^/p{l - C)y 

As usual, the energy and momentum balance equations 
do not close the problem, and we need an additional re- 
lation between the objects of the theory. For Newtonian 
fluids it is know that in the log-layer W^y and K are 
proportional to each other 

W,y = clK, c^«0.5. (32a) 

For the problem of drag reduction by polymers this ratio 
is also some constant Cp ~ 0.25 (in the maximum drag 
reduction regime). For the bubbly flow, we define in 
the same manner: 

W:,y = clK. (32b) 

Clearly, ^im = and for small C (noninteracting 

bubbles) — oc C. It was reported in [21, 22] that 
is slightly smaller than its Newtonian counterpart; we 
therefore write 

cl=cl{l-pC), (32c) 

with a positive coefficient f3 of the order of imity. We are 
not aware of direct measurements of this form in bub- 
bly flows, but it appears natural to assume that the pa- 
rameter (3 is y-independent in the turbulent log-law re- 
gion. We note that the Cauchy- Schwartz inequality can 
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be used to prove that Wxy < K, meaning that all the 



ratios c < 1. 

N.B.P — 



V. DRAG REDUCTION IN BUBBLY FLOWS 

In this section we argue that bubble oscillations are 
crucial in enhancing the effect of drag reduction. This 
conclusion is in line with the experimental observation 
of [12] where bubbles and glass spheres were used under 
similar experimental conditions. Evidently, bubble de- 
formations can lead to the compressibility of the bubbly 
mixture. This is in agreement with the simulation of [13] 
where a strong correlation between compressibility and 
drag reduction were found. 

To make the point clear we start with the analysis of 
the energy balance equation (31). The additional stress 
tensor has a diagonal and an off diagonal part. The 
off-diagonal part has a viscous part that is negligible for 
high Re. The other term can be evaluated using the 
estimate (29), leading to the contribution 

{^{w-U){w-U) -.Vu) ^2{uu:Vu) . (33) 

The expression on the RHS is nothing but the spatial tur- 
bulent energy flux which is known to be very small in the 
log-layer compared to the production term on the RHS of 
Eq. (31). We will therefore neglect the off-diagonal part 
of the stress tensor in the energy equation. The analysis 
of the diagonal part of the stress depends on the issue 
of bubble oscillations and we therefore discuss separately 
oscillating bubbles and rigid spheres. 



A. Drag reduction with rigid spheres 

Consider first situations in which R — 0. This is the 
case for bubbles at small We, or when the bubbles are re- 
placed by some particles which are less dense than neat 
fluid [12]. When the volume of the bubbles is fixed, the 
incompressibility condition for the Newtonian fluid is un- 
changed, and Sii = 0. The diagonal part of r, due to 
the incompressibility condition su = 0, has no contri- 
bution to (TijSij). The energy balance equation is then 
unchanged compared to the Newtonian fluid. The mo- 
mentum balance equation is nevertheless affected by the 
bubbles. Putting (32) into (31), we have 



xy 



p{l - C) 



(34) 



To assess the amount of drag reduction we will consider 
an experiment [14] in which the velocity profile (and thus 
S) is maintained fixed. Drag reduction is then measured 
by the reduction in the momentum flux P. We then have 



P = 



p{l-C + aC)cl 



(35) 



where k is the von-Karman constant. If there are no bub- 
bles (C = 0), the Newtonian momentum flux P^^ reads 



P 

The percentage of drag reduction can be defined as 



(36) 



%DR 



{l-C + aC)cl 



P c" 



(1 - a -h 3/3) C . 



(37) 



Here we assumed that /? ^ 1. At small Re, a = and 
the amount of drag reduction increases linearly with C . 
If Re is very large we expect a « 2, and then the drag is 
enhanced. This result is in pleasing agreement with the 
experimental data in [12]. Indeed, the addition of glass 
beads with density less than water caused drag reduction 
when Re is small, whereas at Re ~ (10^), the drag was 
slightly enhanced. 



B. Drag reduction with flexible bubbles 

If the value of We is sufficiently large such that R ^ 
0, the velocity field is no longer divergenceless. To see 
how this affects the energy cqiiation we consider a single 
bubble with volume V. From the continuity equation: 



/ 



u-dA = V . 



(38) 



Assumed that the bubble is small enough such that the 
velocity field does not change much on the scale of R, 
then we have 



(39) 



V R 

Therefore, the last term in (21) can be approximated as 



4:11— = —V • U 

^ R 3 



(40) 



Next we substitute Eq. (21) into Eq. (24). For smafl 
amplitude oscillations we can neglect the terms propor- 
tional to R^ [23]. The expression for the stress tensor 
simplifies to 



27 



(41) 

For large We, the term p{wi — Ui){wj — Uj)/2 becomes 
larger than the terms Pg — 2'y/R + p. Using Eq. (39) 



P 



4 1 

-^^lS^j6ij + -{Wi - Ui){Wj - Uj) 



(42) 



The extra turbulent dissipation due to the bubble is 
(TijSij). In light of the smallness of the term in Eq. (33) 
we find 



(njs^j) = i-psl) . 



(43) 
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Note that A is an unknown parameter that should de- 
pend on We, and so its value is different in different ex- 
periments. The percentage of drag reduction for various 
values of A are shown in Fig. 1 where we chose a = 2 and 
for simplicity we estimate f3 = 0. One sees that for a = 2 
and A = (where the latter is associated with rigid bub- 
bles), we only find drag enhancement. For small value of 
A, or small amplitudes of oscillations, small concentra- 
tions of bubbles lead (for a = 2) to drag enhancement, 
but upon increasing the concentration we find modest 
drag reduction. Larger values of A lead to considerably 
large degrees of drag reduction. For A = 0.15, the re- 
sult agrees reasonably with Legner's model which pre- 
dicts %DR « 1 - 5(1- C)2/4 [11]. Note that according 
to Legner, there should be considerable drag enhance- 
ment when C = 0. This is of course a nonsensical result 
that is absent in our theory. For ^ = 0.8, %DR « 4C for 
small C. This is the best fit to the experimental results 
which are reported in [14]. 
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FIG. 1: Predicted values of drag reduction with a = 2 and 

different values of A. In a dashed line we reproduce the predic- 
tions of Legner's model which suffer from an unpfiysical drag 
enliancemeut at C = 0, For ^ = (rigid spheres) we find only 
increasing drag enhancement as a function of C. For small 
values of A we have first a slight drag enhancement, and then 
modest drag reduction. For large values of A, associated with 
strong bubble oscillations, we find significant values of drag 
reduction. 



The term |/xs^j is of the same form as the usual dissipa- 
tion term jiSijSij and therefore we write this as: 

|3/2 j^3/2 



4 2x_j^PM 



A 



y 



Vp(1 - C)3/2y ' 



(44) 



where A is an empirical constant. Finally, the energy 
equation becomes 



6(1 -C) + AC K^/'^ 

7^(r^cp7^^ 



W^yS 



(45) 



As before, we specialize the situation to an experiment 
in which S is constant, and compute the momentum flux 



p(l-C)-{l-C + aC) cl 



The degree of drag reduction is then 



24 \ 



(46) 



(47) 



VI. SUMMARY AND DISCUSSION 

The main conclusion of this study is that bubble os- 
cillations can contribute decisively to drag reduction by 
bubbles in turbulent flows. In agreement with the experi- 
mental findings of [12], we find that rigid bubbles tend to 
drag enhance, and the introduction of oscillations whose 
amplitude is measured by the parameter A (Fig.l) in- 
creases the efficacy of drag reduction. 

It is also important to recognize that bubble oscilla- 
tions go hand- in-hand with the compressibility V-it 7^ 0. 
In this sense we are in agreement with the proposition 
of [13] that drag reduction by bubbles is caused by the 
compressibility. There is a difference, however: in [13] 
the flow is free (having only one wall) whereas in our 
case we have a channel in mind. The mechanism of [13] 
cannot appear in our case. On the other hand [13] does 
not allow for bubble oscillations. The bottom line is that 
in both cases the bubble dynamics leads to the existence 
of compressibility, and the latter contributes to the drag 
reduction. 

One drawback of the present study is that the bub- 
ble concentration is taken uniform in the flow. In real- 
ity a profile of bubble concentration may lead to even 
stronger drag reduction if placed correctly with respect 
to the wall. A consistent study of this possibility calls for 
the consideration of buoyancy and the self-consistent so- 
lution of the bubble concentration profile. Such an effort 
is beyond the scope of this paper and must await future 
progress. 

Finally, it should be noted that we neglected the effects 
of viscosity in equations (30) and (31) as we assumed the 
value of Re to be large. For moderately large Re, one 
can take the viscosity effects into account as suggested in 
[18]. 
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